Временные ряды. Основные понятия.Примеры

Будем говорить, что временной ряд это последовательность случайных величин (одномерных или многомерных) \(x_t=(x_{t_1},...,x_{t_n})\), наблюдаемая в моменты времени \(t_1,...,t_n\). Другими словами временной ряд \(x_t\) -это случайный процесс с дискретным временем. Если интервалы времени \(t_k-t_{k-1}=const\) для всех \(k\), то индекс \(t_k\) заменим просто на порядковый номер. Итак случайный процесс \(x_t=(x_1,...,x_n)\) с дискретным временем будем называть временным рядом. В теории обычно предполагается, что временной ряд имеет бесконечную длину и динамика его рассматривается от \(-\infty\) до \(\infty\). В этом случае это записывают так

\({x_t:t = 0,\pm1,\pm2,....}\)

Временные ряды часто изображают графически, по оси абцисс располагают время, наблюдаемые значения процесса отображают по оси ординат. В библиотеке TSA содержатся многие функции небходимые для анализа временных рядов, там же находятся многочисленные примеры временных рядов. Например, временной ряд количество осадков в дюймах за год в Лос Анжелесе

library("TSA")
## Loading required package: leaps
## Loading required package: locfit
## locfit 1.5-9.1    2013-03-22
## Loading required package: mgcv
## Loading required package: nlme
## This is mgcv 1.8-9. For overview type 'help("mgcv-package")'.
## Loading required package: tseries
## 
## Attaching package: 'TSA'
## The following objects are masked from 'package:stats':
## 
##     acf, arima
## The following object is masked from 'package:utils':
## 
##     tar
data(larain) 
plot(larain,ylab='Inches',xlab='Year',type='o')

Другой пример это ряд количества пассажиро-милей в США по месячно

data(airmiles)
plot(airmiles)

Последний ряд интересен тем, что на нем хорошо заметны компоненты из которых состоит ряд. Предполагается, что временной ряд состоит из следующих компонент

\(x_t=m_t+s_t+i_t+u_t+\epsilon_t\)

где через \(m_t\) обозначают тренд - основную тенденцию в динамике ряда, в рассматриваемом примере основная тенденция -это рост. Через \(s_t\) обозначен сезонный эффект, или просто сезонность. Это некоторый эффект характерный для динамики ряда и повторяющийся,вообще говоря, с известным периодом. В данном примере период сезонности равен 12 - число месяцев в году. Через \(i_t\) будем обозначать интервенции- резкое изменение в динамике процесса, вызванное часто внешним воздействием. В примере это известный всем террористический акт 11 сентября 2001 года
Последнюю компоненту ряда , \(u_t+\epsilon_t\), будем называть стационарным случайным процессом. Стогое определение этого понятия будет дано в курсе позднее. А также станет понятно почему оно состоит из двух
частей \(u_t\) и \(\epsilon_t\). Последняя носит название “белый шум”. Почему шум имеет белый цвет и какого другого цвета бывают шумы также позднее будет пояснено. Увидеть компоненту \(u_t+\epsilon_t\) в ряде пока трудно, мешают предшествующие. Деление на компоненты \(m_t,s_t,i_t и u_t+\epsilon_t\) весьма условно, в ходе анализа сезонность или интервенция сами вдруг станут трендом. Модель

\(x_t=m_t+s_t+i_t+u_t+\epsilon_t\)

называют аддитивной моделью.

Возможно представления ряда в виде произведения соответствующих компонент

\(x_t=m_t*s_t*i_t*(u_t+\epsilon_t)\)

В этом случае она носит название мультипликативной модели, которая сводится к аддитивной путем формального логарифмирования.

Простейшие описательные статистики

Для случайного процесса временого ряда \({x_t:t = 0,\pm1,\pm2,....}\) опаределим математическое ожидание.

\(E[x_t]=\mu_t,t=0,\pm1,\pm2,...\)

В общем случае \(\mu_t\) различно для каждого \(t = 0,\pm1,\pm2,....\)

Дисперсия ряда

\(D[x_t]=E(x_t-\mu_t)^2= \sigma_t^2 ,t=0,\pm1,\pm2,...\)

Очень важную роль будет играть автоковариационная функция \(c(x_t,x_s)\) и автокорреляционная функция \(\gamma(x_t,x_s)\) Определяются они следующим образом соответсвенно

\(c(x_t,x_s) =cov(x_t,x_s)=E(x_t-\mu_t)(x_s-\mu_s):s,t=0,\pm1,\pm2,...\)

и

\(\gamma(x_t,x_s) =\frac{c(x_t,x_s)}{\sqrt{ D[x_t]}\sqrt {D[x_s]}}:s,t=0,\pm1,\pm2,...\)

По выборке исторических данных \(x_1,...,x_n\) оценка среднего, дисперсии и автоковариационной функции осуществляется по формулам:

математическое ожидание

\(E[x_t]\approx\overline{x}= \frac{1}{n}\sum_{t=1}^nx_t\)

дисперсия

\(D[x_t]\approx \frac{1}{n-1}\sum_{t=1}^n(x_t- \overline x)^2\)

автоковариационная функция

\(с(k=t-s)\approx \frac{1}{n-k}\sum_{t=1}^{n-k}(x_t- \overline x)(x_{t+k}- \overline x)\)

Для выборки

s<- 3
sigma  <- 2
x <- rnorm(100,s,sigma)
matplot(x, type ="l",col = "blue",lwd = 2)

mean(x)
## [1] 2.902229
var(x)
## [1] 5.065893
acf(x,col="blue",lwd=2)

Белый шум и случайное блуждание

Пусть \(\epsilon_1,...,\epsilon_n\) последовательность независимых нормально распределенных случйаных величин, или так называемый нормальный “белый шум” с математическим ожиданием \(E[\epsilon_t]=0\) и дисперсией \(D[\epsilon_t]=\sigma_{\epsilon}^2:t=1,..n\)

sigma <- 1
e <- rnorm(100,0,sigma)
matplot(e,type = "l",main= "Simulated white noise",lwd = 3,col = "blue", xlab= "time",ylab = "white noise")

Построим процесс \(x_t\) следующим образом

\(x_1=\epsilon_1\)

\(x_2=x_1+\epsilon_1\)

\(....\)

\(x_n=x_{n-1}+\epsilon_n\)

Процесс \(x_t\) называется случайным блужданием. Нетрудно убедиться, что

\(E[x_t]=E[x_{t_1}+\epsilon_t]=E[\epsilon_1+...+\epsilon_t]= 0\),

а

\(D[x_t]=D[x_{t_1}+\epsilon_t]=D[\epsilon_1+...+\epsilon_t]=t\sigma_{\epsilon}^2\).

Пусть \(1\le t\le s\) .Автокорреляционная функция, так как \(E[\epsilon_t=0]\)

\(c(x_t,x_s)=E[(\epsilon_1+...\epsilon_t)(\epsilon_1+...\epsilon_s)]=\) \(\sum_{i=1}^t\sum_{j=1}^sE[\epsilon_i\epsilon_j]\)

Так как \(E[\epsilon_i\epsilon_j]=\sigma_\epsilon^2\) при \(i=j\) и равно 0 при \(i\ne j\) тогда

\(c(x_t,x_s)=t\sigma_\epsilon^2\)

Более того легко видеть,что при \(1\le t\le s\)

\(c(x_t,x_s)=c(x_s,x_t) = t\sigma_\epsilon^2\)

Аналогично легко может быть вычислена и автокорреляционная функция для процесса \(x_t\) случайного блуждания

\(\rho(x_t,x_s)=\frac{c(x_t,x_s)}{\sqrt{D[x_t]D[x_s]}}=\sqrt{\frac{t}{s}}\)

Случайное блуждание на рисунке выглядит следующим образом

x <- diffinv(e,lag = 1,differences = 1,0)
matplot(x,type = "l",main= "Simulated random walk",lwd = 3,col = "blue",xlab ="time",ylab = "random walk")

Скользящее среднее

Пусть \(\epsilon_1,...,\epsilon_n\) по-прежнему процесс “белого шума” Построим процессс \(x_t:t=1,...,n\) по следующему правилу

\(x_t(2)=\frac{\epsilon_t+\epsilon_{t-1}}{2}\)

Простейшие вычисления нам дают, что при всех \(t=1,...n\) \(\mu_t=E[x_t]= E[\frac{\epsilon_t+\epsilon_{t-1}}{2}]= 0\)

a

\(\sigma_{x,t}^2=D[x_t]= D[\frac{\epsilon_t+\epsilon_{t-1}}{2}]=\frac{D[\epsilon_t]+D[\epsilon_{t-1}]}{4}]= 0.5 \sigma_{\epsilon}^2\)

Также нетрудно убедиться,что

\(cov(x_t,x_{t-1})= 0.25 \sigma_{\epsilon}^2\)

а

\(cov(x_t,x_{t-2})= 0\)

Объединяя последние 3 выражения, получим, что автоковариационная функция процесса \(x_t\) равна

\(с(x_t,x_s)=0.5\sigma_{\epsilon}^2\) при \(|t-s|= 0\)

\(с(x_t,x_s)=0.25\sigma_{\epsilon}^2\) при \(|t-s|=1\)

\(с(x_t,x_s)= 0\) при \(|t-s|>1\)

для автокорреляционной функции роцесса \(x_t\)

\(\rho(x_t,x_s)=1\) при \(|t-s|= 0\)

\(\rho(x_t,x_s)=0.5\) при \(|t-s|=1\)

\(\rho(x_t,x_s)= 0\) при \(|t-s|>1\)

Аналогично можно ввести скользящее среднее для трех,четырех и любого произвольного числа \(m\) слагаемых

\(x_t(m)=\frac{\epsilon_t+\epsilon_{t-1}+...+\epsilon_{t-m}}{m}\)

Вычислить математическое ожидание, дисперсию, автоковариационную и автокорреляционную функции процесса \(x_t(m)\) Для любого \(m = 2,3,4,...\)

\(\mu_t=E[x_t]= E[\frac{\sum_{i=0}^{m-1}\epsilon_{t-i}}{m}]= 0\)

\(\sigma_t^2=D[x_t]= D[\frac{\sum_{i=0}^{m-1}\epsilon_{t-i}}{m}]=\frac{\sum_{i=0}^{m-1}D[\epsilon_{t-i}]}{m^2}=\frac{m\sigma_{\epsilon}^2}{m^2}=\frac{\sigma_{\epsilon}^2}{m}\)

Автокоррелляционная функция

\(\rho(x_t,x_s)=1\) при \(|t-s|= 0\)

\(\rho(x_t,x_s)\ne0\) при \(|t-s|<m\)

\(\rho(x_t,x_s)= 0\) при \(|t-s|\ge m\)

x_2<- filter(e,rep(1/2,2))
x_10<- filter(e,rep(1/10,10))

matplot(cbind(e,x_2,x_10),type = "l",lty=1,main= "Simulated white noise and moving 
        averages",lwd = 3,col = c("blue","red","green"),xlab ="time",ylab = "moving averages")
legend("topright", c("White noise","SMA(2)","SMA(10)"),col = c("blue","red","green"))

##Стационарность##

Важной характеристикой временного ряда является свойство стационарности

Определение 1. Cтрогая стационарность.

Случайный процесс c дискретным временем \(x_t:t =0,\pm1,\pm2,...\) называется строго стационарным, если для любого целого \(k>0\)для любых моментов времени \(t_1,....,t_k\) для любого действительного \(s\) распределение процесса в моменты времени \(t_1,....,t_k\) и \(t_1+s,....,t_k+s\) остается неизменным. Формально запишем это как

\(f(x_{t_1},....,x_{t_k}=f(x_{t_1+s},....,x_{t_k+s})\)

Для \(k=1\) из строгой стационарности процесса \(x_t\) следует, что

\(E[x_t]=\mu=const,t=0,\pm1,\pm2,...\)

Дисперсия также постоянна во времени

\(D[x_t]=E(x_t-\mu_t)^2= \sigma^2=const ,t=0,\pm1,\pm2,...\)

При \(k=2\) из строгой стационарности процесса \(x_t\) следует, что автоковариационная функция и автокорреляционная функция для любых \(t\le s\) зависят только от сдвига \(|s-t|\)

\(c(x_t,x_s)=c_{|s-t|}\)

\(\gamma(x_t,x_s)=\gamma_{|s-t|}\)

Определение 2. Стационарность в широком смысле.

Случайный процесс c дискретным временем \(x_t:t =0,\pm1,\pm2,...\) называется стационарным в широком смысле , если математическое ожидание и дисперсия процесса постоянны во времени \(E[x_t]=\mu=const,t=0,\pm1,\pm2,...\)

\(D[x_t]=E(x_t-\mu_t)^2= \sigma^2=const ,t=0,\pm1,\pm2,...\)

Рассмотренные выше примеры процессов белого шума и скользящего среднего для белого шума это примеры строго стационарных процессов

Импорт экономических временных рядов

На сайте http://www.finam.ru/ в разделах “Про рынок->Экспорт данных” и “Про рынок -> Календарь статистики” находится большая коллекция постоянно обновляющихся временных рядов. Которые можно импортировать через csv файл. На языке R программа импорта дневных минутных данных цен на акции Газпром выглядит так

gazprom <- read.csv("c:/tmp/gazprom.csv",header = TRUE)
gaspromClose <- gazprom$X.CLOSE.
gazpromTime <- gazprom$X.TIME.
summary(gazprom)
##  X.TICKER.      X.PER.      X.DATE.        X.TIME.       X.OPEN.     
##  GAZP:522   Min.   :1   02/05/16:522   10:01:00:  1   Min.   :133.9  
##             1st Qu.:1                  10:02:00:  1   1st Qu.:134.5  
##             Median :1                  10:03:00:  1   Median :135.1  
##             Mean   :1                  10:04:00:  1   Mean   :135.1  
##             3rd Qu.:1                  10:05:00:  1   3rd Qu.:135.6  
##             Max.   :1                  10:06:00:  1   Max.   :136.1  
##                                        (Other) :516                  
##     X.HIGH.          X.LOW.         X.CLOSE.         X.VOL.      
##  Min.   :134.0   Min.   :133.9   Min.   :133.9   Min.   :    40  
##  1st Qu.:134.5   1st Qu.:134.4   1st Qu.:134.5   1st Qu.: 11125  
##  Median :135.2   Median :135.1   Median :135.1   Median : 25820  
##  Mean   :135.1   Mean   :135.0   Mean   :135.1   Mean   : 46540  
##  3rd Qu.:135.7   3rd Qu.:135.6   3rd Qu.:135.6   3rd Qu.: 53333  
##  Max.   :136.2   Max.   :136.1   Max.   :136.1   Max.   :639940  
## 
plot(gazpromTime,gaspromClose,type = "b",main = " Gazprom (1 Minute)",lwd = 2,lty = 1)

Импорт статистических данных по торговому балансу России (ежемесячно). Временные ряды удобно при выдаче на график хранить в виде объекта типа time series. Внимание!!! После импорта обязательно построить график и убедиться визуально, что в данных нет явных ошибок

trbalans <- read.csv("c:/tmp/balans.csv",header = TRUE)
summary(trbalans)
##          Äàòà          Ïåðèîä         Ôàêò        Ïðîãíîç       
##  01.04.2005:  1   ßíâ ' 03:  2   Min.   : 0.018   Mode:logical  
##  02.02.2004:  1   Àâã ' 00:  1   1st Qu.: 4.516   NA's:227      
##  02.02.2015:  1   Àâã ' 01:  1   Median : 9.539                 
##  04.11.2003:  1   Àâã ' 02:  1   Mean   : 9.358                 
##  04.12.2003:  1   Àâã ' 03:  1   3rd Qu.:13.787                 
##  05.03.2007:  1   Àâã ' 04:  1   Max.   :20.533                 
##  (Other)   :221   (Other) :220                                  
##       Fact            Åä..èçì.  
##  Min.   : 0.000   ìëðä USD:132  
##  1st Qu.: 4.478   ìëðä.$  : 95  
##  Median : 9.539                 
##  Mean   : 9.465                 
##  3rd Qu.:13.787                 
##  Max.   :52.846                 
## 
trbalansFact <- trbalans$Fact

trbalansTS <- ts(trbalansFact, start = c(1997, 1), frequency = 12)

plot(trbalansTS,type = "b",main = " Trade balans(monthly)",col = "blue",lwd = 2,lty = 1)